Phonon spectrum and bonding properties of Bi 2 Se 3 : Role of strong spin-orbit 

interaction 
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Phonon dispersions of one typical three-dimensional topological insulator Bi2Se3 have been stud- 
ied within density functional theory. The soft modes of two acoustic branches along the Z—F 
and T—F directions within the pure local density approximation will transit to show imaginary 
frequency oscillating after including the spin-orbit interaction (SOI). Similar phenomenon has also 
been observed for Bi2Te3. Besides, we have found that the weak van der Waals forces between two 
Sel layers in Bi2Se3 are strengthened by turning on the SOI. 
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Topological insulators (TIs) have an energy gap in 
their bulk but, due to SOI, possess one or more robust 
metallic states on their edge or surface protected by time- 
reversal symmetry. |j| This novel quantum state, which 
can show various topological quantum effects combined 
with its potential applications involving future spintronic 
devices, quantum computing, and photonics, has at- 
tracted great attention in condensed-matter physics and 
materials science recently. After revealing the topo- 
logical nature of the surface states in the first three- 
dimensional (3D) TI: Bii-^Sba; by angle-resolved pho- 
toemission spectroscopy (ARPES), the second generation 
materials: Bi2Se3, Bi2Te3, and Sb2Te3 were theoretically 
predicted to be the simplest 3D TIs whose surface states 
consist of a single Dirac cone at the Gamma point @. 
The topological nature of these material has already been 
established by recent experiments. [__-__[ 

Here we focus our sight on Bi 2 Sc;j, in which an ARPES 
study [_} discovered a single surface electron pocket with 
a Dirac point below the Fermi level. Actually, before the 
renewing interests due to its topological surface states, 
many efforts have been concentrated on this material ow- 
ing to its large thermoelectric effect. Early X-ray diffrac- 
tion experiment [l(| demonstrated that Bi 2 Se3 possesses 
layered structure along the hexagonal z axis [see Fig. 
1(a)] with alternating layers of Bi and Se in Sel-Bi-Sc2- 
Bi-Sel order (known as quintuple layers). Transport and 
optical experiments [lll - ll3l ] have determined the insult- 
ing band gap to be approximately 0.25-0.35 eV. Raman 
and far-infrared [l_| as well as inelastic neutron scat- 
tering (l5| experiments have been conducted to investi- 
gate the phonon spectrum and phonon density of state 
(DOS). Ultrafast carrier and phonon dynamics tem- 
perature dependent optical characters in a magnetic field, 
[13| Faraday effect [131, an d Raman response [18| exper- 
iments of Bi2Se3 have been studied in recent days. The 
easy cleavage between its neighboring quintuple layers is 
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FIG. 1: (a) Hexagonal crystal structure of Bi2Se3 where a 
quintuple layer with Sel-Bi-Se2-Bi-Sel sequence is presented 
by the red square, (b) Trigonal primitive cell of Bi2Se3. 



due to the weak van der Waals forces between two Sel 
layers. Strong covalent bonding character within each 
quintuple layer has already been analyzed by performing 
density functional theory (DFT) calculations. [ll|] Our 
chemical bonding analysis will show that by including 
SOI the Sel-Sel bonds can be viewed more strengthened 
than that described by pure DFT. The dynamical unsta- 
ble feature of Bi 2 Sc3 and Bi 2 Te 3 is observed by investi- 
gating the phonon dispersions. 

In this letter, the first-principles DFT calculations on 
the basis of the frozen-co re p rojected augmented wave 
(PAW) method of Blochl [2fj] are performed within the 
Vienna ab initio simulation package (VASP)I21(, where 
the local density approximation (LDA) [22| and the 
Perdew, Burke, and Ernzerhof (PBE) [23| form of the 
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TABLE I: Optimized structural parameters and insulating 
band gap (E g ) for Bi2Se3 without and with SOI. For compar- 
ison, experimental values are also listed. 



Method 


a 


c 


u 


V 


E 9 




(A) 


(A) 


(A) 


(A) 


(eV) 


LDA 


4.102 


27.865 


0.4015 


0.2090 


0.0 


LDA+SOI 


4.112 


27.611 


0.4022 


0.2074 


0.3 


Expt. 


4. 143 ° 


28.636° 


0.4008° 


0.2117° 


0.25-0.35° 



° Reference [10|. 6 References [liT - floT ] 




FIG. 2: Band structure for Bi2Se3 without (a) and with (b) 
SOI. The Fermi energy level is set at zero. 
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FIG. 3: Calculated phonon dispersion curves, phonon DOS, 
phonon free energy (F p h), entropy (S), and specific heat at 
constant volume (CV) for Bi2Se3 without and with SOI. Ex- 
perimental results of phonon DOS from Ref [l5[ and Cv from 
Ref [2S| (hollow circles) are also presented. 
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generalized gradient approximation (GGA) are employed 
to describe electron exchange and correlation. To obtain 
accurate total energy, a cutoff energy of 300 eV is used 
for the plane- wave set with a 10x10x10 k point- meshes 
in the Brillouin zone (BZ) for the rhombohedral crystal 
structure. All atoms arc fully relaxed until the Hellmann- 
Feynman forces becoming less than 0.001 eV/A. The Bi 
5<i 10 6s 2 6p 3 and the Se 4s 2 4p 4 orbitals are included as va- 
lence electrons. SOI is included to compare with calcu- 
lations without considering SOI. Within LDA+SOI for- 
malism, the symmetry of the system is switched off. 

Bi2Se3 crystallizes in rhombohedral crystal structure 
with space group i?3m (No. 166) with Bi in 2c(/i, /i, 
At), Sel in 2c(u, u, v), and Se2 in la(0, 0, 0) Wyckoff 
positions. Using GGA and GGA+SOI, we obtain non- 
rational overestimations of the equilibrium volume and 
the band gap. But LDA and LDA+SOI can give proper 
results of structural parameters in comparison with ex- 
periments (see Table I) and also can successfully pre- 
dict Bi 2 Se 3 to be a topological insulator, as shown in 
Fig. 2. Our calculated insulating band gap of about 0.3 
eV within LDA+SOI formalism agrees well with recent 
infrared reflectance and transmission measurements [l3| 
and previous calculations @, [l9| . 

Phonon spectrum has tight relation with dynamical 
stability, phase transition, thermoelectric effect, and su- 
perconductivity. In the present work, we compute the 
phonon spectrum for Bi2Se3 with and without consider- 
ing SOI. Phonon frequency calculations are performed by 
using the supercell approach within the FROPHO code. 
24| To reach high accuracy, we use 3x3x3 rhombohedral 
supercell containing 135 atoms and 2x2x2 Monkhorst- 
Pack fc-point mesh for the BZ integration. The calcu- 
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FIG. 4: Calculated phonon dispersion curves and phonon 
DOS for Bi2Te3 without and with SOI. Experimental phonon 
DOS is taken from Ref [H. 



lated phonon curves along T—Z—F—T—L directions and 
the phonon density of states (DOS) are displayed in Fig. 
3. For comparison, experimental phonon DOS from Ref. 
is also presented. 
In the overall view of Fig. 3, most phonon curves at op- 
tical branches calculated by LDA+SOI are a little lower 
than that without including SOI. Two lowest acoustic 
branches exhibit soft oscillating modes along the Z— F 
and T—F directions within LDA. After including the SOI 
effect, this kind of soft modes transit to show imaginary 
phonon frequencies, which illustrates dynamic unstabil- 
ity of Bi2Se3 single crystal at low temperature. This 
unstable nature has not been reported before, although 
some experimental works [Til - filll lilj ] have focused on its 
lattice vibration. To exclude the specificality of our re- 
sults about Bi2Se3, we have also calculated the phonon 
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TABLE II: Phonon frequencies at the T point for Bi2Se3 and Bi2Te3 in unit of THz. 
results and experimental values are also listed. 



For comparison, previous calculation 



Bi 2 Se 3 



Bi 2 Te 3 



Modes symmetry 



LDA 


LDA+SOI 


LDA+SOI" 


Expt." 


LDA 


LDA+SOI 


GGA+sor 


Expt. c 


1.32 


1.24 


1.24 


1.17 


1.28 


1.27 


1.06 


1.03 


2.26 


2.32 


2.26 


2.20 


1.96 


1.99 


1.62 


1.86 


4.31 


4.15 


4.11 


3.98 


2.38 


3.20 


2.88 


3.05 


5.45 


5.27 


5.13 


5.26 


4.26 


4.04 


3.81 


4.02 


2.64 


2.41 


2.41 


2.04 


2.14 


1.75 


1.45 




4.06 


3.92 


3.91 


3.75 


3.03 


2.87 


2.74 




4.31 


4.16 


4.12 


3.87 


3.15 


2.82 


2.85 




4.93 


4.83 


4.84 


4.80 


3.85 


3.62 


3.56 
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curves and phonon DOS of Bi2Te3 (see Fig. 4). Clearly, 
the imaginary phonon frequencies along T—F direction 
after switching on SOI can be seen. The good agree- 
ment between our calculated phonon DOSs and the ex- 
perimental results further guarantees the correctness of 
our calculations. We notice that after structural analysis 
a recent high-pressure work of Bi2Te3 [25j | owed the a- 
Bi2Te3— >/3-Bi2Te3 phase transition mechanism to a par- 
allel movement of BiTe6 along the [110] direction. Our 
phonon results support their conclusion and supply di- 
rect vibration images for both Bi2Se3 and Bi2Te3. The 
similarity of the phonon results for these two typical 3D 
TIs suggests same phase transition path upon compres- 
sion. In Table II, we present the eigenfrequencies of A\ g , 
AT g , E\, E l g \ E l u , E 1 *, Al u , and Af u phonon modes at 
Gamma point for Bi2Se3 and Bi2Tc3 within LDA and 
LDA+SOI approaches. Good agreement with previous 
calculation (2(| and experiments [l3|, EH EH EH H3] can 
be found. Compared with LDA, LDA+SOI improves the 
agreement with experiments. As for the phonon DOS, 
relative strength of three major peaks and two positions 
of them are successfully reproduced in comparison with 
experiment. The second major peak is relatively higher 
than experiment. These differences between our calcu- 
lations and experiments can be partially attributed to 
the temperature effects. After all, the experiments were 
conducted under 300 and 77 K [H, EBJ. Besides, recent 
experiment (l6l | indicated that the air exposure can signif- 
icantly affect the carrier and phonon dynamics in Bi2Se3 
crystals. The temperature-induced electron-phonon in- 
teractions and defect-induced charge trapping should be 
responsible for its dynamic stability. 

Based on phonon DOS, we further obtain some ther- 
modynamic properties using the quasiharmonic approx- 
imation (28j . under which the Helmholtz free energy 
F( T, V) at temperature T and volume V can be ex- 
pressed as 

F(T, V) = E{V) + F p h (T, V) + F el (T, V), (1) 

where E( V) is the ground-state total energy, F p h(T, V) is 
the phonon free energy and F e i(T, V) is the thermal elec- 
tronic contribution. For present material, we focus only 




FIG. 5: Valence charge density of Bi2Se3 in the az-plane with- 
out (a) and with SOI (b). Contour lines are drawn from 0.0 
to 0.1 at 0.01 e/au 3 intervals. 



on the contribution of atom vibrations. The F p h(T,V) 
can be calculated by 



F ph (T,V) = k B T s (w)ln 



2 sinh 



2k B T 



dw, (2) 



where H and k B are Planck and Boltzmann constants, re- 
spectively, uj represents the phonon frequencies and g(uS) 
is the phonon DOS. Here, we only present results at equi- 
librium volumes from both LDA and LDA+SOI. The en- 
tropy (S) can be determined by S=-T(d F/d T) v . The 
specific heat at constant volume Cy can be directly cal- 
culated through 



( Huj 



exp 



k B T 



\k B Tj (exp^-1) 



-du. (3) 



Calculated results of F p h(T), S, and Cy and previous 
experimental values of Cy are shown in Fig. 3, where it 
can be found that the calculated Cy agrees well with ex- 
periment [29| . No evident differences have been observed 
between LDA and LDA+SOI calculations. We hope that 
our results can provide instruction for further study. 

To understand the chemical-bonding characters of 
Bi2Se3, we investigate the valence charge density distri- 
bution by plotting plane charge density in the iz-plane 
in both LDA and LDA+SOI cases (see Fig. 5) as well 
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TABLE III: Calculated charge and volumes according to Bader partitioning as well as the bond lengthes for Bi2Se3 without 
and with SOI. 



Method 


Qs(Bi) 


Q B (Sel) 


Q B (Se2) 


V s (Bi) 


V s (Sel) 


V B (Se2) 


Sel-Bi 


Se2-Bi 


Sel-Sel 




(e) 


(e) 


(e) 


(A 3 ) 


(A 3 ) 


(A 3 ) 


(A) 


(A) 


(A) 


LDA 


14.05 


6.55 


6.79 


24.38 


29.16 


28.25 


2.84 


3.04 


3.35 


LDA+SOI 


14.07 


6.54 


6.78 


24.61 


28.35 


28.18 


2.84 


3.04 


3.27 



as by performing the Bader analysis [3(| (see Table III). 
The Sel— Bi, Se2— Bi, and Sel— Sel bond distances are 
also tabulated in Table III. The correlated minimum val- 
ues of charge density along the Sel— Bi, Se2— Bi, and 
Sel— Sel bonds are presented in Fig. 5. Obviously, our 
present LDA results demonstrate that the binding be- 
tween Bi and Se layers exhibits predominant covalency 
and the interaction between two Sel layers is of the Van 
der Waals type, which coincides well with previous con- 
clusions, [la Hi Comparing the results obtained by LDA 
and LDA+SOI, we find that the main influence on the 
structure, electronic distribution, and also the bonding 
nature of SOI lies between two Sel layers. The weak 
van der Waals forces bonding of Sel-Sel is strengthened 
by including SOI. The z axis of the hexagonal crystal 
and the bond length of Sel-Sel is shortened by SOI (see 
Tables I and III). As shown in Table III, the atomic vol- 
ume of Sel is diminished largely by LDA+SOI in com- 
parison with the pure LDA calculations, while the in- 
crease or decrease amplitudes of Vs(Bi) and Vs(Se2) 
are not so evident. With SOI, the ionic charges of Bi 2 Sc 3 
also have a little change compared with LDA. While Bi 
atom losses 0.95 electrons within LDA, it transfer 0.93 
electrons to Se within LDA+SOI formalism. The ionic 
charges can be represented as Bi 2 K95+ Sel 2 l - 55 ~Se2 - 79 - 



and Bi°- 93+ Sel^ 54 -Se2 o ra - within LDA and LDA+SOI, 
respectively. In two approaches, each Se2 atom gains 
more electrons than Sel atom. 

In summary, we have demonstrated that the anti- 
crossing feature in the band structure around the T point 
of Bi 2 Se3 can be successfully predicted by full geome- 
try optimization within LDA and LDA+SOI formalisms. 
Good agreement of structural parameters and the eigen- 
frequencies of phonon modes at T point for Bi2Se3 and 
Bi2Te3 between our theoretical results and experiments 
can be found. Within LDA+SOI, we have observed imag- 
inary phonon frequencies along Z—F and T—F direc- 
tions for Bi 2 Se3 and along T~F direction for Bi 2 Te 3 , 
which supports recent high- pressure work [25| . Their dy- 
namic stability in experiments can be regarded as due 
to the temperature-induced electron-phonon interactions 
and defect-induced charge trapping. We have found that 
the SOI plays an important role in describing the weak 
van der Waals forces between two Sel layers. 
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dations for Development of Science and Technology of 
China Academy of Engineering Physics under Grant No. 
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